{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Check that gradients are computed correctly \n",
    "(compare numeric derivative approximation computed using official C implementation and derivative computed by backpropagation using Pytorch version)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "import numpy as np\n",
    "import pandas as pd\n",
    "import subprocess\n",
    "\n",
    "import matplotlib.pyplot as plt\n",
    "plt.rcParams['figure.dpi']=200\n",
    "\n",
    "import torch\n",
    "import torch.nn.functional as F\n",
    "\n",
    "from vmaf_torch import VMAF, yuv_to_tensor, tensor_to_yuv\n",
    "\n",
    "np.set_printoptions(precision=5, suppress=True)\n",
    "torch.set_printoptions(precision=5, sci_mode=False)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "'cuda'"
      ]
     },
     "execution_count": 2,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "device = 'cuda' if torch.cuda.is_available() else 'cpu'\n",
    "device"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "torch.Size([1, 1, 1080, 1920])"
      ]
     },
     "execution_count": 3,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# load reference yuv as tensors\n",
    "#yuv_path_ref = \"/storage/data/NFLX_dataset_public/ref/BigBuckBunny_25fps.yuv\"\n",
    "yuv_path_ref = \"/storage/data/NFLX_dataset_public/ref/BirdsInCage_30fps.yuv\"\n",
    "\n",
    "width = 1920\n",
    "height = 1080\n",
    "num_frames = 1  \n",
    "\n",
    "ref_y, ref_u, ref_v = yuv_to_tensor(yuv_path_ref, width, height, num_frames, color='yuv')\n",
    "ref_y = ref_y.to(device)\n",
    "ref_u = ref_u.to(device)\n",
    "ref_v = ref_v.to(device)\n",
    "ref_y.shape"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [],
   "source": [
    "# initialize VMAF Pytorch version\n",
    "vmaf = VMAF(NEG=False, enable_motion=True, clip_score=True)\n",
    "vmaf = vmaf.to(device)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# initialize VMAF C version\n",
    "class VMAF_C():\n",
    "    '''Utility class for calling reference vmaf executable'''\n",
    "\n",
    "    def __init__(self, vmaf_executable=\"vmaf\", vmaf_model_version=\"default\", verbose=True):\n",
    "        self.vmaf_executable = vmaf_executable\n",
    "        self.verbose = verbose\n",
    "\n",
    "        if vmaf_model_version == \"default\":\n",
    "            self.vmaf_model_param = \"\"\n",
    "        elif vmaf_model_version == \"NEG\":\n",
    "            self.vmaf_model_param = \"--model version=vmaf_v0.6.1neg\"\n",
    "\n",
    "    def table_from_path(self, ref_path, dist_path, width, height, num_frames):\n",
    "        vmaf_out_csv_path = 'vmaf_out.csv'\n",
    "        vmaf_param = f\"{self.vmaf_executable} -r {ref_path} -d {dist_path} -w {width} -h {height} --frame_cnt {num_frames} -p 420 -b 8 --threads 16 -q --csv -o {vmaf_out_csv_path} {self.vmaf_model_param}\"\n",
    "        if self.verbose:\n",
    "            print('Executing:', vmaf_param)\n",
    "        p = subprocess.run(vmaf_param.split(), stdout=subprocess.PIPE, stderr=subprocess.PIPE)\n",
    "        if self.verbose:\n",
    "            print('Reading:', vmaf_out_csv_path)\n",
    "        df = pd.read_csv(vmaf_out_csv_path)\n",
    "\n",
    "        return df\n",
    "\n",
    "    def score_from_path(self, ref_path, dist_path, width, height, num_frames):\n",
    "        df = self.table_from_path(ref_path, dist_path, width, height, num_frames)\n",
    "        score = df['vmaf'].iloc[:num_frames].mean()\n",
    "\n",
    "        return score\n",
    "\n",
    "    def table_from_tensors(self, ref_tup, dist_tup):\n",
    "        # TODO check shapes and rounding here\n",
    "        ref_save_path = './ref.yuv'\n",
    "        dist_save_path = './dist.yuv'\n",
    "        width = ref_tup[0].shape[-1]\n",
    "        height = ref_tup[0].shape[-2]\n",
    "        num_frames = ref_tup[0].shape[0]\n",
    "        if self.verbose:\n",
    "            print('Saving tensors to:', ref_save_path, 'and', dist_save_path)\n",
    "        tensor_to_yuv(*ref_tup, yuv_path=ref_save_path)\n",
    "        tensor_to_yuv(*dist_tup, yuv_path=dist_save_path)\n",
    "\n",
    "        return self.table_from_path(ref_save_path, dist_save_path, width, height, num_frames)\n",
    "\n",
    "    def score_from_tensors(self, ref_tup, dist_tup):\n",
    "        df = self.table_from_tensors(ref_tup, dist_tup)\n",
    "        score = df['vmaf'].mean()\n",
    "\n",
    "        return score\n",
    "\n",
    "\n",
    "vmaf_c = VMAF_C(vmaf_executable='vmaf', vmaf_model_version='default')  # wrapper class for vmaf executable"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "tensor([[[[0.04000, 0.04000, 0.04000, 0.04000, 0.04000],\n",
       "          [0.04000, 0.04000, 0.04000, 0.04000, 0.04000],\n",
       "          [0.04000, 0.04000, 0.04000, 0.04000, 0.04000],\n",
       "          [0.04000, 0.04000, 0.04000, 0.04000, 0.04000],\n",
       "          [0.04000, 0.04000, 0.04000, 0.04000, 0.04000]]]], device='cuda:0',\n",
       "       requires_grad=True)"
      ]
     },
     "execution_count": 6,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "#kernel_size = 3\n",
    "kernel_size = 5\n",
    "\n",
    "weight = torch.ones((1, 1, kernel_size, kernel_size))/kernel_size**2\n",
    "\n",
    "weight = weight.to(device)\n",
    "weight.requires_grad = True\n",
    "\n",
    "weight"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [],
   "source": [
    "ref_y_conv = F.conv2d(ref_y, weight=weight, padding='same')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Max pixel value: 230.2799835205078 min pixel value: 19.959999084472656\n"
     ]
    }
   ],
   "source": [
    "with torch.no_grad():\n",
    "    print('Max pixel value:', ref_y_conv.max().item(), 'min pixel value:', ref_y_conv.min().item())\n",
    "\n",
    "if ref_y_conv.max()>255 or ref_y_conv.min()<0:\n",
    "    print('WARNING: images have pixels out of [0,255] range, the approximate gradients might be incorrect!')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "f = 61.04484939575195\n"
     ]
    }
   ],
   "source": [
    "f = vmaf(ref_y, ref_y_conv)\n",
    "print('f =', f.item())"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "tensor([[[[148.16345, 168.11429, 179.23662, 177.40149, 162.86739],\n",
       "          [166.76881, 188.79459, 199.59097, 194.09003, 174.89992],\n",
       "          [176.37427, 198.75858, 207.17685, 199.56595, 177.43134],\n",
       "          [174.69728, 194.16537, 200.52112, 190.49844, 168.83699],\n",
       "          [163.28206, 178.01871, 180.62872, 170.31906, 150.94447]]]],\n",
       "       device='cuda:0')"
      ]
     },
     "execution_count": 10,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# compute derivatives via backpropagation\n",
    "f.backward()\n",
    "weight.grad"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Derivative computed via backpropagation using Pytorch:  207.17684936523438\n"
     ]
    }
   ],
   "source": [
    "# select one element of the filter\n",
    "coef_address = (0,0,kernel_size//2,kernel_size//2)\n",
    "#coef_address = (0,0,1,2)\n",
    "#coef_address = (0,0,0,0)\n",
    "\n",
    "derivative = weight.grad[coef_address].item()\n",
    "\n",
    "print(\"Derivative computed via backpropagation using Pytorch: \", derivative)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(tensor([[[[0.04000, 0.04000, 0.04000, 0.04000, 0.04000],\n",
       "           [0.04000, 0.04000, 0.04000, 0.04000, 0.04000],\n",
       "           [0.04000, 0.04000, 0.05000, 0.04000, 0.04000],\n",
       "           [0.04000, 0.04000, 0.04000, 0.04000, 0.04000],\n",
       "           [0.04000, 0.04000, 0.04000, 0.04000, 0.04000]]]], device='cuda:0'),\n",
       " tensor([[[[0.04000, 0.04000, 0.04000, 0.04000, 0.04000],\n",
       "           [0.04000, 0.04000, 0.04000, 0.04000, 0.04000],\n",
       "           [0.04000, 0.04000, 0.03000, 0.04000, 0.04000],\n",
       "           [0.04000, 0.04000, 0.04000, 0.04000, 0.04000],\n",
       "           [0.04000, 0.04000, 0.04000, 0.04000, 0.04000]]]], device='cuda:0'))"
      ]
     },
     "execution_count": 12,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "eps = 0.01\n",
    "\n",
    "weight_plus_eps = weight.detach().clone()\n",
    "weight_plus_eps[coef_address] += eps \n",
    "\n",
    "weight_minus_eps = weight.detach().clone()\n",
    "weight_minus_eps[coef_address] -= eps \n",
    "\n",
    "weight_plus_eps, weight_minus_eps"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "metadata": {},
   "outputs": [],
   "source": [
    "ref_y_conv_plus_eps = F.conv2d(ref_y, weight=weight_plus_eps, padding='same')\n",
    "ref_y_conv_minus_eps = F.conv2d(ref_y, weight=weight_minus_eps, padding='same')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 14,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Max pixel value: 232.58998107910156 min pixel value: 20.299999237060547\n"
     ]
    }
   ],
   "source": [
    "print('Max pixel value:', ref_y_conv_plus_eps.max().item(), 'min pixel value:', ref_y_conv_plus_eps.min().item())\n",
    "if ref_y_conv_plus_eps.max()>255 or ref_y_conv_plus_eps.min()<0:\n",
    "    print('WARNING: images have pixels out of [0,255] range, the approximate gradients might be incorrect!')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 15,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Max pixel value: 227.96998596191406 min pixel value: 19.619998931884766\n"
     ]
    }
   ],
   "source": [
    "print('Max pixel value:', ref_y_conv_minus_eps.max().item(), 'min pixel value:', ref_y_conv_minus_eps.min().item())\n",
    "if ref_y_conv_minus_eps.max()>255 or ref_y_conv_minus_eps.min()<0:\n",
    "    print('WARNING: images have pixels out of [0,255] range, the approximate gradients might be incorrect!')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 16,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Saving tensors to: ./ref.yuv and ./dist.yuv\n",
      "Executing: vmaf -r ./ref.yuv -d ./dist.yuv -w 1920 -h 1080 --frame_cnt 1 -p 420 -b 8 --threads 16 -q --csv -o vmaf_out.csv \n",
      "Reading: vmaf_out.csv\n",
      "Saving tensors to: ./ref.yuv and ./dist.yuv\n",
      "Executing: vmaf -r ./ref.yuv -d ./dist.yuv -w 1920 -h 1080 --frame_cnt 1 -p 420 -b 8 --threads 16 -q --csv -o vmaf_out.csv \n",
      "Reading: vmaf_out.csv\n",
      "f_plus=63.114191 f_minus=58.959366\n",
      "Derivative computed numerically using C implementation:  207.74124999999975\n"
     ]
    }
   ],
   "source": [
    "# compute approximate derivative using reference C implementation\n",
    "f_plus = vmaf_c.score_from_tensors((ref_y, ref_u, ref_v), (torch.round(ref_y_conv_plus_eps), ref_u, ref_v))\n",
    "f_minus = vmaf_c.score_from_tensors((ref_y, ref_u, ref_v), (torch.round(ref_y_conv_minus_eps), ref_u, ref_v))\n",
    "derivative_approx = (f_plus - f_minus) / (2 * eps)\n",
    "print(f'{f_plus=}', f'{f_minus=}')\n",
    "\n",
    "print(\"Derivative computed numerically using C implementation: \", derivative_approx)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 17,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Derivative computed numerically using C implementation:  207.74124999999975\n",
      "Derivative computed via backpropagation using Pytorch:  207.17684936523438\n",
      "Difference -0.5644006347653772\n"
     ]
    }
   ],
   "source": [
    "print(\"Derivative computed numerically using C implementation: \", derivative_approx)\n",
    "print(\"Derivative computed via backpropagation using Pytorch: \", derivative)\n",
    "print(\"Difference\", derivative-derivative_approx)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 18,
   "metadata": {},
   "outputs": [],
   "source": [
    "# sanity check: compute approximate derivative using Pytorch VMAF implementation\n",
    "#f_plus = vmaf(ref_y, ref_y_conv_plus_eps)\n",
    "#f_minus = vmaf(ref_y, ref_y_conv_minus_eps)\n",
    "#grad_approx = (f_plus - f_minus) / (2 * eps)\n",
    "#print(f'{f_plus=}', f'{f_minus=}')\n",
    "#print(f'{grad_approx=}')\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 19,
   "metadata": {},
   "outputs": [],
   "source": [
    "# sanity check: compute approximate derivative using Pytorch VMAF implementation after rounding\n",
    "#f_plus = vmaf(ref_y, torch.round(ref_y_conv_plus_eps))\n",
    "#f_minus = vmaf(ref_y, torch.round(ref_y_conv_minus_eps))\n",
    "#grad_approx = (f_plus - f_minus) / (2 * eps)\n",
    "#print(f'{f_plus=}', f'{f_minus=}')\n",
    "#print(f'{grad_approx=}')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 20,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Coefficient: 0 0\n"
     ]
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Derivative computed numerically using C implementation:  147.5766\n",
      "Derivative computed via backpropagation using Pytorch:  148.1634521484375\n",
      "Difference 0.5868521484374867\n",
      "Coefficient: 0 1\n",
      "Derivative computed numerically using C implementation:  168.33855\n",
      "Derivative computed via backpropagation using Pytorch:  168.11428833007812\n",
      "Difference -0.2242616699218729\n",
      "Coefficient: 0 2\n",
      "Derivative computed numerically using C implementation:  180.22220000000004\n",
      "Derivative computed via backpropagation using Pytorch:  179.2366180419922\n",
      "Difference -0.985581958007856\n",
      "Coefficient: 0 3\n",
      "Derivative computed numerically using C implementation:  177.47409999999988\n",
      "Derivative computed via backpropagation using Pytorch:  177.4014892578125\n",
      "Difference -0.07261074218737917\n",
      "Coefficient: 0 4\n",
      "Derivative computed numerically using C implementation:  162.4544499999999\n",
      "Derivative computed via backpropagation using Pytorch:  162.8673858642578\n",
      "Difference 0.41293586425791773\n",
      "Coefficient: 1 0\n",
      "Derivative computed numerically using C implementation:  166.35649999999984\n",
      "Derivative computed via backpropagation using Pytorch:  166.76881408691406\n",
      "Difference 0.4123140869142219\n",
      "Coefficient: 1 1\n",
      "Derivative computed numerically using C implementation:  189.91985000000005\n",
      "Derivative computed via backpropagation using Pytorch:  188.79458618164062\n",
      "Difference -1.1252638183594286\n",
      "Coefficient: 1 2\n",
      "Derivative computed numerically using C implementation:  199.11280000000033\n",
      "Derivative computed via backpropagation using Pytorch:  199.59097290039062\n",
      "Difference 0.47817290039029103\n",
      "Coefficient: 1 3\n",
      "Derivative computed numerically using C implementation:  194.31735000000003\n",
      "Derivative computed via backpropagation using Pytorch:  194.09002685546875\n",
      "Difference -0.2273231445312831\n",
      "Coefficient: 1 4\n",
      "Derivative computed numerically using C implementation:  175.58760000000007\n",
      "Derivative computed via backpropagation using Pytorch:  174.89991760253906\n",
      "Difference -0.6876823974610033\n",
      "Coefficient: 2 0\n",
      "Derivative computed numerically using C implementation:  176.16064999999992\n",
      "Derivative computed via backpropagation using Pytorch:  176.374267578125\n",
      "Difference 0.2136175781250813\n",
      "Coefficient: 2 1\n",
      "Derivative computed numerically using C implementation:  198.33694999999983\n",
      "Derivative computed via backpropagation using Pytorch:  198.75857543945312\n",
      "Difference 0.4216254394532939\n",
      "Coefficient: 2 2\n",
      "Derivative computed numerically using C implementation:  207.74124999999975\n",
      "Derivative computed via backpropagation using Pytorch:  207.17684936523438\n",
      "Difference -0.5644006347653772\n",
      "Coefficient: 2 3\n",
      "Derivative computed numerically using C implementation:  199.17909999999992\n",
      "Derivative computed via backpropagation using Pytorch:  199.56594848632812\n",
      "Difference 0.3868484863282049\n",
      "Coefficient: 2 4\n",
      "Derivative computed numerically using C implementation:  177.87240000000003\n",
      "Derivative computed via backpropagation using Pytorch:  177.43133544921875\n",
      "Difference -0.4410645507812774\n",
      "Coefficient: 3 0\n",
      "Derivative computed numerically using C implementation:  174.83085000000003\n",
      "Derivative computed via backpropagation using Pytorch:  174.69728088378906\n",
      "Difference -0.133569116210964\n",
      "Coefficient: 3 1\n",
      "Derivative computed numerically using C implementation:  192.54184999999993\n",
      "Derivative computed via backpropagation using Pytorch:  194.16537475585938\n",
      "Difference 1.6235247558594494\n",
      "Coefficient: 3 2\n",
      "Derivative computed numerically using C implementation:  200.5752000000001\n",
      "Derivative computed via backpropagation using Pytorch:  200.5211181640625\n",
      "Difference -0.05408183593760896\n",
      "Coefficient: 3 3\n",
      "Derivative computed numerically using C implementation:  190.76625000000007\n",
      "Derivative computed via backpropagation using Pytorch:  190.49844360351562\n",
      "Difference -0.2678063964844455\n",
      "Coefficient: 3 4\n",
      "Derivative computed numerically using C implementation:  167.55964999999975\n",
      "Derivative computed via backpropagation using Pytorch:  168.8369903564453\n",
      "Difference 1.2773403564455634\n",
      "Coefficient: 4 0\n",
      "Derivative computed numerically using C implementation:  162.39854999999997\n",
      "Derivative computed via backpropagation using Pytorch:  163.2820587158203\n",
      "Difference 0.8835087158203407\n",
      "Coefficient: 4 1\n",
      "Derivative computed numerically using C implementation:  176.8408000000001\n",
      "Derivative computed via backpropagation using Pytorch:  178.01870727539062\n",
      "Difference 1.1779072753905382\n",
      "Coefficient: 4 2\n",
      "Derivative computed numerically using C implementation:  180.7120000000001\n",
      "Derivative computed via backpropagation using Pytorch:  180.62872314453125\n",
      "Difference -0.08327685546885277\n",
      "Coefficient: 4 3\n",
      "Derivative computed numerically using C implementation:  170.19565000000014\n",
      "Derivative computed via backpropagation using Pytorch:  170.31906127929688\n",
      "Difference 0.12341127929673235\n",
      "Coefficient: 4 4\n",
      "Derivative computed numerically using C implementation:  149.46184999999977\n",
      "Derivative computed via backpropagation using Pytorch:  150.94447326660156\n",
      "Difference 1.4826232666017916\n"
     ]
    }
   ],
   "source": [
    "# do all of this for every coefficient\n",
    "\n",
    "vmaf_c.silent = True \n",
    "\n",
    "derivatives_c = []\n",
    "derivatives_pytorch = []\n",
    "dif = []\n",
    "\n",
    "for i in range(kernel_size):\n",
    "    for j in range(kernel_size):\n",
    "        coef_address = (0,0,i,j)\n",
    "        print('Coefficient:', i, j)\n",
    "        \n",
    "        derivative = weight.grad[coef_address].item()\n",
    "\n",
    "        weight_plus_eps = weight.detach().clone()\n",
    "        weight_plus_eps[coef_address] += eps \n",
    "        weight_minus_eps = weight.detach().clone()\n",
    "        weight_minus_eps[coef_address] -= eps \n",
    "\n",
    "        ref_y_conv_plus_eps = F.conv2d(ref_y, weight=weight_plus_eps, padding='same')\n",
    "        ref_y_conv_minus_eps = F.conv2d(ref_y, weight=weight_minus_eps, padding='same')\n",
    "\n",
    "        if ref_y_conv_plus_eps.max()>255 or ref_y_conv_plus_eps.min()<0:\n",
    "            print('WARNING: images have pixels out of [0,255] range, the approximate gradients might be incorrect!')\n",
    "        if ref_y_conv_minus_eps.max()>255 or ref_y_conv_minus_eps.min()<0:\n",
    "            print('WARNING: images have pixels out of [0,255] range, the approximate gradients might be incorrect!')\n",
    "            \n",
    "        # compute approximate derivative using reference C implementation\n",
    "        f_plus = vmaf_c.score_from_tensors((ref_y, ref_u, ref_v), (torch.round(ref_y_conv_plus_eps), ref_u, ref_v))\n",
    "        f_minus = vmaf_c.score_from_tensors((ref_y, ref_u, ref_v), (torch.round(ref_y_conv_minus_eps), ref_u, ref_v))\n",
    "        derivative_approx = (f_plus - f_minus) / (2 * eps)\n",
    "\n",
    "        print(\"Derivative computed numerically using C implementation: \", derivative_approx)\n",
    "        print(\"Derivative computed via backpropagation using Pytorch: \", derivative)\n",
    "        print(\"Difference\", derivative-derivative_approx)\n",
    "        \n",
    "        derivatives_c.append(derivative_approx)\n",
    "        derivatives_pytorch.append(derivative)\n",
    "        dif.append(derivative-derivative_approx)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 21,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(0.5739042109375305, 0.4530605698159361)"
      ]
     },
     "execution_count": 21,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "dif = np.abs(np.array(dif))\n",
    "dif.mean(), dif.std()"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "th",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.10.5"
  },
  "orig_nbformat": 4
 },
 "nbformat": 4,
 "nbformat_minor": 2
}
